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Abstract. Gamma-ray burst (GRB) afterglows are well described by synchrotron emission originating from the interaction 
between a relativistic blast wave and the external medium surrounding the GRB progenitor. We introduce a code to reconstruct 
spectra and light curves from arbitrary fluid configurations, making it especially suited to study the effects of fluid flows 
beyond those that can be described using analytical approximations. As a check and first application of our code we use it 
to fit the scaling coefficients of theoretical models of afterglow spectra. We extend earlier results of other authors to general 
circumburst density profiles. We rederive the physical parameters of GRB 970508 and compare with other authors. 

We also show the light curves resulting from a relativistic blast wave encountering a wind termination shock. From high 
resolution calculations we find that the observed transition from a stellar wind type light curve to an interstellar medium type 
light curve is smooth and without short-time transitory features. 
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INTRODUCTION 

In the fireball model, Gamma-Ray Burst (GRB) afterglows are thought to be the result of synchrotron radiation 
generated by electrons during the interaction of a strongly collimated relativistic jet from a compact source with 
its environment (for recent reviews, see l22l[T8ll ). Initially the resulting spectra and light curves have been modelled 
using only the shock front of a spherical explosion and a simple power law approximation for the synchrotron radiation 
(e.g. ll27llT7ll24ll ). One or more spectral and temporal breaks were used to connect regimes with different power law 
slopes. For the dynamics the self similar Blandford-McKee (BM) approximation of a relativistic explosion was used 
Oj. These models have been refined continuously. More details of the shock structure were included (e.g. ll6l [Toll ). 
more accurate formulae for the synchrotron radiation were used (e.g. ll26ll ) and efforts have been made to irnplement 
collimation using various analytical approximations to the jet structure and lateral spreading behavior (see \9] for an 
overview). On to p of that, there have been studies focussing on arrival time effects (e.g. |12]) and some numerical 
simulations (e.g. lilllH). 

In this paper we introduce a code to reconstruct spectra and light curves from AMRVAC, a high performance 
relativistic hydrodynamics code ifTsll . We verify our method by applying it to the analytically well understood BM 
solution. Because different authors have recently started using the circumstellar density structure as a fitting parameter 
when fitting the BM solution to afterglow data|25], we generalize existing scaling coefficient prescriptions from the 
literature ig] from insterstellar medium (ISM, for which the inverse radial slope of the density distribution k is zero) 
and stellar wind (k — 2) to general k. These scaling coefficients are tabulated in the appendix and can be directly 
used when fitting to afterglow data. We finish this part of the paper by comparing fit results to GRB 970508 using our 
prescription to those of other authors. 

Following this, we apply our radiation code to study the visible effect of the blast wave encountering a wind 
termination shock. Our simulations, done at high resolution to make sure we accurately probe the timescales at which 
the encounter is expected to take place, confirm the prediction of ifTgll of a smooth transition between two power law 
regimes in the observerd light curve. 

With the exception of the wind termination shock section, most of the work presented in this paper is also presented 
inH. 



DESCRIPTION OF THE RADIATION CODE 



The code takes as input a series of snapshots of relativistic hydrodynamics configurations on a grid. The grids represent 
a spherically symmetric fluid configuration and all grid cells are assumed to emit a fraction of their energy as radiation. 



This fraction of course has to be small enough not to affect the dynamics, since the post-processing approach does not 
allow for feedback. For the time being we restrict ourselves to the optically thin case. In this section and the next we 
will use BM solution for adiabatic expansion of the blast wave to provide the content of the grid snapshots. The BM 
solution takes two input parameters: £52, the explosion energy in units of 10^^ erg and hq, the circumburst number 
density at a characteristic distance of 10^^ cm. 

Four ignorance parameters are provided to the code at runtime: p, ^n, £e and Eb, denoting respectively the slope 
of the relativistic particle distribution, the fraction of particles accelerated to this relativistic distribution at any given 
time, the fraction of thermal energy that is carried by the relativistic electrons and the fraction of thermal energy that 
resides in the (tangled-up) magnetic field. 

In this work we consider synchrotron radiation only. All grid cells contain a macro scop ic number of radiating 
particles and the radiation from these particle distributions is calculated following (l^ and f 23'], but with two important 
differences: the transition to the lab frame is postponed as long as possible and no assumption about the dynamics of 
the system is used anywhere as this should be provided by the snapshot files. 

For the emitted power per unit frequency of a typical electron we have 
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Here qe denotes the electron charge, nie the electron mass and B' the local magnetic field strength. Comoving 
quantities are primed. The function =S contains the local shape of the spectrum, including the local effect of electron 
cooling, resulting in J3°^ v'/^ for values below Vcr.m, in v(i-p)/2 between Va-.m and Vcr,M and an exponential drop 
beyond Vcr,M- incorporates an integration over all pitch angles between electron velocities and the local magnetic 
field and an integration over the accelerated particle distribution. We use a power law particle distribution with a lower 
cut-off Lorentz factor 7,,, and an upper cut-off Lorentz factor Ym- Directly behind the shock front the lower cut-off 
is determined by the downstream density, thermal energy and the ignorance parameters, while the upper cut-off is 
initially set to infinity. The evolution of the accelerated particle distribution when a shocked fluid parcel moves further 
downstream is determined by adiabatic cooling and synchrotron radiation losses. The characteristic frequencies V^^^ 
and V^^„ are related to the bounding Lorentz factors Ym and 7,„ via v' °^ B y^- More details on equation ([T]), the critical 
frequencies and the full shape of Q can be found in ||3|] . 

We emphasize that =S describes the local synchrotron spectrum for a location where all electrons have undergone 
exactly the same amount of cooling. The observed spectrum from the entire blast wave is a superposition of many such 
spectra. The result of this superposition is that the exponential drop beyond the locally differing Vc-.m gets smoothened 
into a steepening of the power law slope of the spectrum by a factor — 1 /2, and the different Vcr,M values together 
determine the position of the cooling break, beyond which this steepening of the slope sets in. 

Assuming isotropic radiation in the comoving frame, we arrive at 

dP-P' 1 dP' 

( v') = — ( v') (2) 

per solid angle 

To get to the received power per unit volume in the lab frame, we have to apply the correct beaming factors, Doppler 
shift the frequency and multiply the above result for a single particle with the lab frame particle density: 

with ji now denoting the cosine of the angle between the fluid velocity and the observer (unprimed, so measured in 
the lab frame), j3 the fluid velocity in units of c and n the number density. 
Finally, the flux the observer receives at a given observer time is given by 

Fiy) = ^[ 4^iv\vmi-[i^)cdAdt,. (4) 
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Here robs is the observer distanc^E approximately the same for all fluid cells (though the differences in arrival times 
are taken into account). The area A denotes the equidistant surface. For every emitting time fg a specific intersecting 
(with the radiating volume) surface exists from which radiation arrives exactly at The integration over the emission 
times te (represented in the different snapshot files) requires an extra beaming factor and a factor of c to transform the 
total integral to a volume integral. 



SCALING COEFFICIENTS AND APPLICATION TO GRB970508 



Especially for high Lorentz factors, the shape of the spectrum is dominated by the radiation coming from a very thin 
slab right behind the shock front. The observed emission from this slab depends on the various model parameters via 
power laws and a heuristic fit function can be constructed with the gradual transition between regimes being handled 
by a free parameter, the sharpness factor s. While the proportionalities of this function can be determined via simple 
scaling arguments, we have to set its scale using our radiation code. Also, in more detailed calculations like those done 
here the gradual transitions are included automatically and we can use this to provide the correct dependence of s on 
p and k. This eliminates s as a free parameter, simplifying the fit to the data and allowing the shape of the transition to 
help determine whether a particular model fits the data or not. Depending on the order of the peak frequency v,„ and 
V,„ we have two options for a fit function. 

If the cooling break lies beyond the peak frequency ('spectrum 1', for easy comparison to |@]), the flux is best 
approximated (i.e. valid up to a few percent) by 
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Here F„, j is the flux at the synchrotron peak frequency for an infinite sharp break (i.e. the point where the two power 
laws describing region D and G would intersect when extrapolated), v,„ \ the synchrotron peak frequency and v^.i the 
cooling frequency. The peak flux and critical frequencies depend on the physical input for the BM solution (explosion 
energy and circumburst density profile) and the ignorance parameters. These dependencies are summarized in the 
appendix. The sharpness functions i and are given by 

i„,,i =2.2-0.52/7, ScA = l.6-Q3%p-Q.l6k + Q.Qlipk. (6) 

When the order of the breaks is reversed ('spectrum 5') the smooth power law for both breaks is given by 
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where Fc^i denotes the peak flux for infinite sharpness 5 and the prescriptions for the sharpness are 

ic,5 =0.66-0. 16/t, i,„,5 =3.7-0.94/:> + 3.64fc-1.16pL 



(7) 



(8) 



Once again these are valid up to a few percent. Given their accuracies, all sharpness prescriptions are consistent with 

a. 

Various authors have used flux scaling equations to derive the physical properties of GRB 970508 from afterglow 
data lUK^H^dll]. This provides us with a context to illustrate the scaling laws presented in this paper. We will use the 
fit parameters obtained from broadband modeling by [11]. They have fit simultaneously in time and frequency while 
keeping as a fitting parameter Because the only model dependencies that have been introduced by this approach are 
the scalings of t and v (and no scaling coefficients), their fit results are still fully consistent with our flux equations. 
Using the cosmology £Im = 0.27, = 0.73 and Hubble parameter Hq — 11 km s^' Mpc^', they have for the 



' For cosmological distances fobs denotes the luminosity distance and redsliift terms ( 1 + z) need to be inserted in tlie appropiate places in the 
equations. 



observer distance in units of 10 cm robs2?< — 1-635 and redshift z = 0.835 [16], leading, at fo^ ^/ — 23.3 days, to 
Vc,i = 9.21 • 10^3 Hz, v„,,i = 4.26 • lO'" Hz, F„i.i = 0.756 mJy, p = 2.22 and k = 0.0307. 

Both flT] and fs] take for the hydrogen mass fraction of the circumburst medium X = 0.7, which in our flux 
equations is mathematically equivalent (though conceptually different) to setting = ( 1 +X) /2 = 0.85. Unfortunately 
this still leaves us with four variables to determine (e^, ££, £52, no) and only three constraints (peak flux, cooling and 
peak frequency). From a theoretical study of the microstructure of collisionless shocks 1 13] arrives at Ee s/eb and 
we include this as an additional constraint to obtain a closed set of equations. 

For the values quoted above we obtain: £52 = 0.155, no = 1.28, £b ~ 0.1057, £e = 0.325. For comparison we give 
some of the values obtained by other authors. |i5J obtain for the ISM case: £52 = 3.5, no — 0.03, Eg — 0.09, Ee =0.12. 
ifj] obtain for the ISM case: £52 = 0.12, no = 22, Eb = 0.012, Ee = 0.57. Both use p = 2.2. Finally [11] obtain for 
k = 0.0307: £52 — 0.435, no ~ 0.0057, Eb = 0.103, Ee = 0.105. This comparison serves to emphasize the sensitivity 
of the inferred blast waves physics on the scaling coefficients of the model. 



WIND TERMINATION SHOCK ENCOUNTER 
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FIGURE 1. Results for wind termination shock encounter. The left figure shows a comoving density snapshot of the shock 
profile during the encounter at emission time tg = 5.49 • 10^ s. The schematic description of |21] is shown for comparison and a 
clear difference between densities further downstream is visible. It should be noted that all different areas of the shock are resolved. 
For example, the forward shock region (smallest, rightmost region) is resolved by ^ 30 cells. The right figure shows the resulting 
light curve at 5 • lO'^ Hz, with 100 data points devoted to 0.3 - 1 day and 100 datapoints to the following 19 days. A smooth 
transition to the power law behaviour corresponding to a BM shock wave expanding into a homogeneous environment is visible. To 
get complete coverage of all observer times, radiation calculated directly from the BM solution for the blast wave at Lorentz 200 
down to 23 has been added to the observed flux. 

Different predictions exist in the literature concerning the observable effect of a relativistic blast wave encounte ring 
a wind-termination shock. A short transitory feature is predicted by ["^Pl, whereas predict a smooth transition. fll\ 
present a detailed analysis of the reverse shock behaviour during the encounter, whereas [19] combine an analysis in 
a planar geometry with simulation results. We perform an RHD simulation using AMRVAC and the initial conditions 
suggested by [/I], with enough temporal and spatial resolution to resolve any short transity features occuring on the 
time scales predicted by [21] (i.e. a few hours in observer time). The output of this simulation is then studied using the 
radiation code presented in this paper and jUtl, for an observer frequency of 5 • lO'^ Hz (optical). 

The resulting light curve shows no short transitory feature but a smooth transition instead, confirming lower 
resolution results from |fl9ll , as can be seen from fig.[T] The discrepancy between our results and ll2lll can be attributed 
partly to the differences in the downstream fluid density profiles between the two studies as illustrated in fig.[T] A more 
detailed argument will be presented elsewhere ^ . 

Following is a short summary of the initial conditions of the simulation. We start with a BM blastwave profile at 
shock Lorentz factor 23, for a shock with £52 = 10.0, no = 3, A: = 2. We position the wind termination shock at a 
radial distance of 1 .6 • lO'** cm. At this point the density increase by a factor 4 and remains fixed further outward. The 
ignorance parameters used in the radiation calculation are p = 2.5, (^at = 1 .0, ££ = 0. 1 and Eb = 0.01 . We have ignored 
electron cooling (i.e. kept Jm at infinity) and differ from l2lll in that we do not take the slight increase in the magnetic 



TABLE 1. Constants setting scale of flux 
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field in the reverse shock region into account. This increase (approximately a factor 1 .2^) is not sufficient to explain 
the reported transitory feature and its omission does not alter our conclusions. 

The grid resolution is determined by the number of base cells (120) and the maximum refinement levels (i.e. the 
number of times a cell can be spUt in two to increase resolution, 15 in our simulation). This grid represents a radial size 
of 6 • IQI^ cm, and therefore the pre-break shock width (A/? at Lorentz factor 23) is resolved at an effective resolution 
of 1200 cells. The temporal resolution is 1.556 • 10^ s (the encounter lasts ^ 1.4 • 10^ s in emission time). 
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APPENDIX: SCALING COEFFICIENTS 

Below we provide the scaling coefficients and model fit functions for general k values. The energy £52 has been 
normalized to units of 10^^ erg, the luminosity distance roi,s_28 to 10^^ cm and the observer time tobs4 is expressed in 
days. The leftmost peak fluxes for infinitely sharp transitions are: 

'"oto,28 

= C,.(^)'.^.ei-nr.Ep.t^.il+z)^rnJy. (10) 



The critical frequencies for the different regimes are given by: 
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The flux functions and critical frequencies still contain the coefficients Co, etc. They are listed below, with Cqq etc. 
denoting purely numerical constants, whose values are tabulated in table [T| 
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